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ABSTRACT 

The standard thin disk model predicts that when the accretion rate is over 
a small fraction of the Eddington rate, which corresponds to L > 0.06LEdd, the 
inner region of the disk is radiation-pressure-dominated and thermally unstable. 
However, observations of the high/soft state of black hole X-ray binaries with 
luminosity well within this regime (0.01L Edd < L < 0.5L Edd ) indicate that the 
disk has very little variability, i.e., quite stable. Recent radiation magnetohydro- 
dynamic simulations of a vertically stratified shearing box have confirmed the 
absence of the thermal instability. In this paper, we revisit the thermal stability 
by linear analysis, taking into account the role of magnetic field in the accretion 
flow. By assuming that the field responses negatively to a positive tempera- 
ture perturbation, we find that the threshold of accretion rate above which the 
disk becomes thermally unstable increases significantly compared with the case 
of not considering the role of magnetic field. This accounts for the stability of 
the observed sources with high luminosities. Our model also presents a possible 
explanation as to why only GRS 1915+105 seems to show thermally unstable 
behavior. This peculiar source holds the highest accretion rate (or luminosity) 
among the known high state sources, which is well above the accretion rate 
threshold of the instability. 
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1. Introduction 

The standard thin disk is a milestone in the development of accretion disk theory. It 
successfully explains many observations such as the "big blue bump" in active galactic nuclei 
(see review in Frank et al. 2002) and the spectrum of high/soft state of black hole X-ray 
binaries (Gierlihski & Done 2004; hereafter GD04). However, some questions still remain 
unsolved (e.g., Koratkar & Blaes 1999). One of them, for example, is how to explain the 
observed hard X-ray emission in black hole sources. The temperature of the thin disk is too 
low to produce the X-ray emission. 

In this paper we will discuss another puzzle, which is its thermal stability. Since the 
discovery of the standard thin accretion disk (e.g., Shakura and Sunyaev 1973), many efforts 
have been made to examine its various stability. When the mass accretion rate is higher 
than a few percent of Eddington rate (Msdd = lO^Edd/c 2 ), which roughly corresponds to 
> 0.06LEdd> the innermost region of the disk will be radiation-pressure-dominated. It was 
found that in this case the disk will be both secularly (Lightman & Eardley 1974) and ther- 
mally (Shakura & Sunyaev 1976; Piran 1978) unstable. Time-dependent global numerical 
calculations (Honma et al. 1992; Szuszkiewicz & Miller 1998; Janiuk et al. 2002; Li et 
al. 2007) found that the local thermal instability will result in the "limit-cycle" behavior. 
However, observations of the high/soft state of black hole X-ray binaries have raised doubt 
to the above prediction (GD04). It has been well established that the accretion flow in this 
state is described by the standard thin disk model (Zdziarski & Gierlihski 2004; Done et al. 
2007). The luminosity of the sources compiled in GD04 ranges from 0.01 to 0.5 ^Edd, with 
the highest one even exceeding the theoretical stable limit by nearly one order of magnitude. 
Thus we should expect to see the limit-cycle behavior. In contrast to this expectation, how- 
ever, observations show little variability which convincingly indicates that they are thermally 
stable. 

Some efforts have been made to solve this puzzle. For example, if a large fraction of 
the dissipated energy is channeled into a corona or outflow/jet the disk would be stable 
(Svensson & Zdziarski 1994). However, the sources adopted in GD04 are all disk-dominated. 
The absence of hard X-ray and radio emission directly rules out the existence of corona or 
jet. Another idea is to assume that the viscous stress is proportional to the gas pressure only, 
rather than the sum of the gas and radiation pressure (Lightman & Eardley 1974; Stella & 
Rosner 1984). As pointed out in GD04, the problem with this idea is that the outcome of this 
modified viscosity is in conflict with the observation that the color temperature correction 
is constant for different accretion rates. Moreover, on the theoretical side, in contrast to the 
above two ideas, numerical simulations have shown that only a very small fraction of the 
dissipated energy is channeled into the corona, which is consistent with the absence of hard 
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X-ray emission in high state of black hole X-ray binaries, and the stress is proportional to 
the sum of gas and radiation pressure (Hirose et al. 2006; Hirose et al. 2009a). Recently, 
the issue of thermal instability has been addressed by the radiation magnetohydrodynamic 
simulations of a vertically stratified shearing box (Hirose et al. 2009a). Their results indicate 
that the radiation-dominated disk is stable over ~ 40 cooling timescales. They explained 
the stability as the lack of correlation between the stress and the pressure within the cooling 
timescale. 

In this paper, we revisit the thermal stability of thin disks by linear analysis. Different 
from previous analyses, we include the role of the magnetic pressure and further assume 
that the magnetic field becomes weaker when the gas temperature increases. Our analysis 
indicates that the onset accretion rate for the instability becomes larger compared with the 
previous result. The organization of our paper is as follows. In §2, we deliver the equations 
and conduct an analysis of the thermal instability. The results are shown in §3. The last 
section (§4) is devoted to a summary and discussion. 



2. Accretion disk model with toroidal magnetic fields 

2.1. Basic equations 

The disk structure is described by the following equations. The equation of vertical 
hydrostatic equilibrium is 

dp , dtp 

where p is the total pressure varying with z, p is the mass density, ip is the gravitational 
potential which is assumed to be of the Newtonian form, i.e., ip(R,z) = —GM/^/R 2 + z 2 , 
where M is the mass of the black hole. Under the thin disk approximation, we have dp/dz ~ 
—P tot /H and dip/dz ~ GMH/R 3 , where H is the scale height, and P tot is the total pressure 
at the midplane. Defining the surface density as E = 2pH, we have 

P -P +P +P _ GMY.H 

1 tot 1 gas ~r J- rad ~T~ J- mag 2,R? ' 

Here the gas pressure P gas = 2/cepT /tor, where &b is the Boltzmann constant, mn is the 
hydrogen mass, and T is the temperature at the midplane. The radiation pressure P ra d = 
aT 4 /3, where a is the radiation constant. Since in the innermost region of the accretion 
flow the magnetic field is dominated by its azimuthal component B^, the magnetic pressure 
P mag = B 2 /8n. 

The energy equation is written as the balance between the viscous heating and the 
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cooling through radiation and advection, i.e., 

Qtis = Qra,d + Qadv- (3) 

The vertically integrated viscous heating rate is 

Qt is = -TryRjl, (4) 



where T r(p is the stress (see eq. 8 below), and f2 = fix with Qk = a/ GM/R 3 is the Keplerian 
angular velocity. The radiative cooling rate is 

n- 32(jT4 ^ 

^rad = —3^- > ( 5 ) 

where r = is the vertical optical depth, /t ~ 0.4 cm 2 g -1 , since the opacity is dominated 
by the electron scattering in the inner region of the disk, and a is the Stefan-Boltzmann 
constant. The advective cooling term is given by (Abramowicz et al. 1995) 

^adv" 2nR2 ■ W 

We set £ = 1.5 in the present work. Our calculations indicate that the role of advection term 
is negligible in term of the stability analysis of the radiation-dominated thin disk. This term 
is included for the discussion of the thermal equilibrium solutions in §3.2, which contain 
optically thick advection-dominated accretion disks (or slim disks; Abramowicz et al. 1988). 

The angular momentum conservation equation is expressed as 

M(Q K R 2 - Zin) = 2nR 2 T rv , (7) 



where l in = yj GMR m is the specific angular momentum at the inner edge of the disk, and 
R m is fixed to be 3R g , with R g = 2GM/c 2 being the Schwarzschild radius. 

It is now widely believed that the MHD turbulence associated with the magneto- 
rotational instability (MRI) is the main mechanism for the angular momentum transfer 
in accretion disks (e.g., Balbus & Hawley 1998). MHD simulations have shown that, the 
stress accounting for the angular momentum transport is dominated by the Maxwell stress, 
rather than the Reynolds stress (Hawley et al. 1995). Besides, we have also learnt from 
the simulations that the ratio of the Maxwell stress to the total pressure nearly maintains 
constant (Hawley & Krolik 2001; Machida et al. 2006; Pessah et al. 2007; Hirose et al. 
2009a). We therefore have 

T rv = 2aP tot H, (8) 
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where a is the viscous parameter. 

To close the set of equations, we need one more constraint, which is to connect the 
magnetic field with other physical quantities. As the thermal instability is mainly concerned 
in this study, here, we relate the response of the magnetic field to the perturbation of the 
scale height in a linearized fashion, i.e., 

^ = -7^ (9) 

where 7 is a free parameter. The integrated form of eq. (9) is $ 7 = B^H 1 = constant. We 
assume that 7 > 0, i.e., the magnetic field will become weaker with an increase of height (or 
temperature). This assumption is supported by the MHD numerical simulation of Machida 
et al. (2006) where they found that the magnetic field becomes stronger when the disk 
shrinks vertically, although we note that their simulation is for a hot accretion flow rather 
than a thin disk. We will discuss the physical meaning of <£> 7 in §4. With the above equations 
we can obtain a thermal equilibrium solution at a certain radius for given M, M, a and $ 7 
(§3.2). 



2.2. Thermal instability analysis 

Since the thermal instability timescale is much shorter than the viscous timescale, the 
surface density E is taken to be constant. Defining /3 gas = P gas /Ptot and Anag = -Pmag/-Ptot> 
we can obtain from eqs. (2), (4)- (6) that 

dmP tot =d\nH 

= fi gas (d In T - din H) + 4(1 - f3 gas - /3 mag )d In T + 2/3 mag d In (10) 

and 

d In Q+ is - d ln(Q rad + Q" v ) = din T rv - 4(1 - / adv )rf In T + f adv (d In M + Id In H), (11) 

where the advection factor / adv is defined as / adv = <5 adv /(<5~ d + <5 adv )- Eqs- (7)-(9) give 
that 

dlnM = dlnT rv = dlnP tot + dlnH, (12) 

and 

dlnB v = --fdlnH. (13) 
Combining eqs. (10) and (13), we have 

dlnP tot = dlnH = 4 -y g as-4An a g dlnT (14) 
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Substituting eqs. (3), (12) and (14) into eq. (11), we finally get 



d(Qtis Qr&d Qa,dv) 

dT 



T _ 2 - 5/3 gas - 4(1 + 7 )/3 mag - 6/ adv + 8/ adv # 



gas 



+ 47)/ ad v/3mag 



£ Qtis 



1 + 



'gas 



+ 2 7 A 



.(15) 



HKlg 



Since the thermal instability condition is [d(Q^ is — — <3 adv )/<9T]£ > 0, and the denom- 
inator of eq. (15) is always positive when 7 > 0, then the thermal instability criterion can 
be taken as 



A = 2 - 5/3 gas - 4(1 + 7 )/3 mag - 6/ adv + 8/ adv /3 gas + (8 + 4 7 )/ adv /3 mag > 0. 



(16) 



For the cases without magnetic field, i.e., /3 mag = 0, the above criterion is reduced to (Gu & 
Lu 2007): 

2 — 5/3 gas — 6/ adv + 8/3 gas / adv > 0. 



3. Results 

3.1. Stability of disks with different magnetic field strength 

We calculate the thermal equilibrium solution of the above equations. The left panel 
of Fig. 1 shows the value of A at it! = 10i? g for various M and /3 mag . The value of 7 is 
selected to be 1 as an illustration. We can see that when /3 mag < 0.24, i.e., the magnetic field 
is relatively weak, there exist two critical mass accretion rates corresponding to A = for 
each given /9 mag . When /3 mag > 0.24, we always have A < 0. That means the disk will be 
thermally stable for any M. The results are similar for the cases of other 7. The right panel 
of Fig. 1 shows the variation of the two critical accretion rates with /3 mag , where the solid and 
dashed lines correspond to the lower and upper critical rates, respectively. The lower critical 
rate is of interest to us since it is relevant to the thin disks. It corresponds to the threshold of 
the accretion rate above which the disk is unstable. The upper one corresponds to slim disks 
with / adv ~ 1/3 (see §3.2). For each 7, the region inside the corresponding parabolic curve 
denotes the thermally unstable solutions, while the region outside corresponds to the stable 
ones. The solid lines in the right panel illustrate that the threshold of accretion rate above 
which the disk is thermally unstable increases with the ratio of magnetic pressure to total 
pressure. Also, we find that if 7 takes a larger value, i.e., the magnetic field reacts more 
intensely according to a perturbation in the scale height, the critical rate increases more 
significantly with /3 mag , and thus a smaller cutoff of the possible /9 mag for thermal instability 
appears. 

We interpret the above result as follows. Suppose there is a small increase in the 
temperature T, the radiative cooling will increase following eq. (5), while the response of 
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viscous heating to temperature is not so simple. As the temperature increases, the radiation 
pressure and further the scale height increase accordingly, which will result in the increase 
of viscous heating (ref. eqs. 4 and 8). The key point is that the magnetic field B v will 
decrease (eq. 9), and this effect appears more important when 7 is relatively large. Therefore, 
compared with the case without magnetic field, the increase of the total pressure and thus the 
scale height will become less significant. In other words, the inclusion of magnetic pressure 
will result in a weaker dependence of viscous heating on temperature, thus the magnetic field 
can conspicuously suppress the thermal instability. 



3.2. Local thermal equilibria 

It is convenient to plot an M ~ £ diagram of local disk solutions to directly show the 
properties of disk stability. For simplicity, we also select the case of 7 = 1 as an instance, 
which is presented in Figure 2. On those curves in this figure, a negative slop indicates 
that the solution is thermally unstable. We can see that the threshold of accretion rate for 
thermal instability, corresponding to the lower turning points of each curve, rises accordingly 
with increasing $1, which is in agreement with Fig. 1. For the dotted and the short dashed 
curves, the thresholds of accretion rate are all ~ M Edd , roughly one order of magnitude 
larger than the case without magnetic field (the long dashed curve), which is supposed to 
explain the observed high luminosity of thermally stable high/soft state of some black hole 
X-ray binaries (GD04). While for the solid curve, as the corresponding $1 is large enough, 
the two turning points both disappear. 

All the curves in Fig. 2 share the same upper branch. That is because in the cases 
of large accretion rates, i.e., > lOM^d, the advective cooling becomes important and the 
radiation pressure turns to be very large. Compared with the large radiation pressure, the 
magnetic pressure can be neglected (the value of P mag / P ra d is of the order of ~ 1CT 2 at 
the upper turning point of each curve except the solid one whereas / a d v approaches 1/3, 
and decreases rapidly as the accretion rate increases). Hence, the effect of magnetic field is 
negligible, and all the curves with different $1 approach the slim disk limit. It is qualitatively 
the same when 7 takes other values. We would like to emphasize that there is no requirement 
that the disk solution would trace the solution curves in Figure 2 when the increase of the 
accretion rate is adequately slow. In fact, the conservation of $ 7 may not hold beyond a 
thermal timescale. 
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4. Summary and discussion 

Previous analysis shows that the radiation-dominated standard thin disk is thermally 
unstable when L > 0.06LEdd- This conclusion, however, is in conflict with observations 
of, e.g., the high state of black hole X-ray binaries. Most previous analyses neglect the 
role of magnetic pressure. In this paper, by taking the magnetic pressure into account, we 
have revisited the thermal stability of the radiation-dominated standard thin disk. We have 
derived a general criterion for the thermal instability (eq. 16). We find that the threshold 
of accretion rate above which the solution becomes unstable increases significantly with the 
increasing ratio of magnetic pressure and total pressure, /3 mag (refer to Fig. 1). If /3 mag is 
large enough, say /3 mag > 0.24 for 7 = 1 or /3 mag > 0.12 for 7 = 3, our model predicts that 
the disk will be stable for any M (again refer to Fig. 1). We do not favor this result since 
the required /3 mag is significantly larger than the typical value of /3 mag ~ 0.1 given by general 
numerical simulations. For /3 mag ~ 0.1, from Fig. 1 we expect that the thermal stability 
of the high state of black hole X-ray binaries with luminosity as high as 0.5^^ can be 
explained if 7 > 30. 

The key assumption in our analysis is that during the thermal perturbation the changes 
of magnetic field and of scale height satisfy a constraint described by eq. (9), namely, the 
magnetic field B v will become weaker with the increase of temperature (or equivalently scale 
height H). The parameter 7 in eq. (9) denotes how strong the response of the magnetic 
field is to a thermal perturbation. Therefore, compared with the case of not including the 
magnetic pressure, the increase of total pressure and further the scale height and viscous 
heating with the temperature become weaker (ref. eqs. 2, 4 & 8). This is the reason why 
the disk tends to be thermally stable. 

Unfortunately we are unclear about the value (or the range) of 7 because of the com- 
plexity of processes such as MRI, dynamo, and reconnection which can strengthen or weaken 
the magnetic field within timescales shorter than the thermal one. One extreme case is that 
all these processes are turned off thus the magnetic field is frozen in the accretion flow. In 
this case, we obviously have 7 = 1, and $1 = B V H represents half of the toroidal magnetic 
flux per unit radius. The conservation of $1 is exactly what has been adopted in the ther- 
mal stability analysis of a magnetically dominated "low-/3" disk (Machida et al. 2006). In 
another work studying this type of disk (Oda et al. 2010), the advection rate of toroidal 
magnetic flux (= 2v r B ip H, where v r = M/27rrE is the radial velocity of the accretion flow) 
is assumed to be constant during a thermal perturbation. This corresponds to 7 = 3, which 



throughout the paper we set a = 0.1. If the value of a is smaller, the required 7 can be smaller to 
explain the observed stable high state. 
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is because from eqs. (2), (7) and (8), we can deduce that v r oc H 2 and therefore B v oc H~ 3 . 
We also notice that the model constructed in Oda et al. (2009) corresponds to 7 = 2. In 
this case, the conservation of $2 = B^H 2 comes from the specific configuration of both the 
magnetic field and the outer boundary condition, which does not hold a simple physical 
meaning. In spite of the uncertainty in the value of 7, we believe that the positive sign of 
7, or equivalently the negative response of B^ to T, should capture the real physics and 
have essential influence on the thermal stability of radiation-dominated thin disks. It will 
be interesting to investigate the value of 7 by detailed numerical simulations. On the other 
hand, it is noteworthy that the reversal of toroidal magnetic field was found by vertically 
stratified shearing box simulations of thin disks (e.g., Brandenburg et al. 1995; Stone et al. 
1996), and most recently by global simulations (O'Neill et al. 2010). This phenomenon was 
reported to occur on ~ 10 orbit timescales, comparable to the thermal timescale. If this phe- 
nomena is coupled with the thermal perturbation process, which is unclear to us, it will be a 
problem for our model, since it will be in conflict with our assumption of the conservation of 
$JE Another concern comes when we try to explain the stability of the radiation-dominated 
simulations presented in Hirose et al. (2009b) based on our model. We find that it requires a 
rather large value of 7 to explain some of their simulations. For example, their model 0519b, 
the most radiation-dominated one, requires 7 > 7. Such a high value feels uncomfortable. 
Of course, it may be too demanding to require a simple one-dimensional analytical model to 
completely explain three-dimensional MHD simulations. 

Our model predicts that if the magnetic pressure is not too strong, e.g., /3 mag < 0.24 
for 7 = 1, the sources with accretion rates higher than the threshold should still manifest 
thermal instability. We note in this context that GRS 1915+105 may be one evidence 
for our prediction. Different from other sources, GRS 1915+105 holds strong variability 
which seems to be well interpreted as the thermal instability of a radiation-dominated thin 
disk (e.g., Belloni et al. 1997; Janiuk et al. 2002). GD04 speculated that the reason for 
the uniqueness of this source might be that it goes to higher luminosity (or equivalently, 
accretion rate) than all other sources in their sample. Actually GRS1915+105 is likely to 
be super-Eddington (Done et al. 2004). This prediction could discriminate our model from 
others such as the one presented in Hirose et al. (2009a). 
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2 In this context we note that our model actually only requires the negative response of the absolute value 
of B v to temperature. But even so, if the field reversal is coupled with thermal perturbation process, since 
the timescales of the two processes are comparable, our assumption of eq. (9) will still be incorrect. 



- 10 - 



grant 2009CB824800, the National Natural Science Foundation of China (grants 10821302, 
10825314, 10833002, and 11073015), and the CAS/SAFEA International Partnership Pro- 
gram for Creative Research Teams. 

REFERENCES 

Abramowicz, M. A., Czerny, B., Lasota, J. P., Szuszkiewicz, E. 1988, ApJ, 332, 646 
Abramowicz, M. A., Chen, X. -M., Kato, S., Lasota, J. P., & Regev, O. 1995, ApJ, 438, L37 
Balbus, S. A., & Hawley, J. F. 1998, Rev. Mod. Phys., 70, 1 

Belloni, T., Mendez, M., King, A. R., van der Klis, M., & van Paradijs, J. 1997, ApJ, 479, 
L145 

Brandenburg, A., Nordlund, A., Stein, R. F., & Torkelsson, U. 1995, ApJ, 446, 741 
Done, C, Gierlihski, M., Kubota, A. 2007, A&A Rev., 15, 1 
Done, C, Wardzihski, C, Gierlihski, M. 2004, MNRAS, 349, 393 

Frank, J., King, A., Raine, D. 2002, Accretion Power in Astrophysics (Cambridge: Cam- 
bridge Univ. Press) 

Gierlihski, M., & Done, C. 2004, MNRAS, 347, 885 

Gu, W. -M., & Lu, J. -F. 2007, ApJ, 660, 541 

Hawley, J. F., Gammie, C. F., & Balbus, S. A. 1995, ApJ, 440, 742 
Hawley, J. F., & Krolik, J. H. 2001, ApJ, 548, 348 
Hirose, S., Krolik, J. H., & Stone, J. M. 2006, ApJ, 640, 901 
Hirose, S., Krolik, J. H., & Blaes, O. 2009a, ApJ, 691, 16 
Hirose, S., Blaes, O., & Krolik, J. H. 2009b, ApJ, 704, 781 
Honma, F., Matsumoto, R., & Kato, S. 1992, PASJ, 44, 529 
Janiuk, A., Czerny, B., & Siemiginowska, A. 2002, ApJ, 576, 908 
Koratkar, A., & Blaes, O. 1999, PASP, 111, 1 



Li, S. -L., Xue, L., & Lu, J. -F. 2007, ApJ, 666, 368 

Lightman, A. P., & Eardley, D. M. 1974, ApJ, 187, LI 

Machida, M., Nakamura, K. E., & Matsumoto, R. 2006, PAS J, 58, 193 

Oda, H., Machida, M., Nakamura, K. E. & Matsumoto, R. 2009, ApJ, 697, 16 

Oda, H., Machida, M., Nakamura, K. E. & Matsumoto, R. 2010, ApJ, 712, 639 

O'Neill, S. M., Reynolds, C. S., Miller, M. C, & Sorathia, K. A. 2010. larXiv: 1009. 18821 

Pessah, M. E., Chan, C. K., & Psaltis, D. 2007, ApJ, 668, L51 

Piran, T. 1978, ApJ, 221, 652 

Shakura, N. L, & Sunyaev, R. A. 1973, A&A, 24, 337 
Shakura, N. L, & Sunyaev, R. A. 1976, MNRAS, 175, 613 
Stella, L., & Rosner, R. 1984, ApJ, 277, 312 

Stone, J. M., Hawley, J. F., Gammie, C. F., & Balbus, S. A. 1996, ApJ, 463, 656 

Svensson, R., & Zdziarski, A. A. 1994, ApJ, 436, 599 

Szuszkiewicz, E., & Miller, J. C. 1998, MNRAS, 298, 888 

Zdziarski, A. A., & Gierlihski, M. 2004, Prog. Theor. Phys. Suppl., 155, 99 



This preprint was prepared with the A AS IATgX macros v5.2. 



- 12 - 




Fig. 1. — The thermal stability of the disk at R — 10R g for various M and /3 mag . M = 10M Q 
and a — 0.1 are adopted. The left panel shows the variation of A with M for given values 
of /9 mag , where 7 is set to be 1. The solid line corresponds to the none magnetic field cases. 
The right panel shows the variation of the lower and upper critical accretion rates (A = 0) 
with /3 mag . Each of the curves corresponds to a specific value for 7. 
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L(g • cm 2 ) 

Fig. 2. — The thermal equilibrium curves of the thin disk at 10i? g for different $1, which is 
defined as $1 = B V H and is in unit of Gs-cm. Other model parameters are M = 10M Q and 
a = 0.1 . 



